Michael Hemmer
Note that this is just a very brief introduction to polynomials. For a quick reference we refer to the Wikipedia or for a more elaborate introduction to any class book on elementary algebra.
A polynomial is either zero, or can be written as the sum of one or more non-zero terms. The number of terms is finite. A term consist of a constant coefficient and a monomial, that is, the product of zero or more variables. Each variable may have an exponent that is a non-negative integer. The exponent on a variable in a term is equal to the degree of that variable in that term. A term with no variables is called a constant term. The degree of a constant term is 0.
For example, -7x3y is a term. The coefficient is -7, the monomial is x3y, comprised of the variables x and y, the degree of x is three, and the degree of y is one. The total degree of the entire term is the sum of the degrees in each variable. In the example above, the degree is 3 + 1 = 4.
A one-variable (univariate) polynomial f of degree n has the following form:
| f = anxn + an-1xn-1 + ... + a2x2 + a1x + a0 |
The coefficient a0 is called the constant coefficient, an is called the leading coefficient. If f is not the zero polynomial the leading coefficient is not zero. The polynomial is called monic if an = 1. In case the coefficient domain of f possess a greatest common divisor (gcd) the content of f is the gcd of all coefficients of f. For instance, the content of 12 x3 + 6 is 6.
A multivariate polynomial is a polynomial in more than one variable. According to the number of variables it is possible to further classify multivariate polynomials as bivariate, trivariate etc. In contrast to univariate polynomials the terms of a multivariate polynomial are not completely ordered by their total degree. However, given a certain order on the variables it is possible to define a lexicographic order on the terms. Given this order the leading coefficient of a multivariate polynomial is defined as the coefficient of the highest term. For instance the leading coefficient of the multivariate polynomial p = 5 x3y + 7xy2 is 7, given that y has an higher order than x.
However, it is also possible to interpret a multivariate polynomial as a univariate polynomial in that variable. For instance the trivariate polynomial
| q = x5 + 7x2y1z2 + 13x1y2z2 ∈ ℤ[x,y,z] |
| q = (13x1y2 + 7x2y1)z2 + x5z0 ∈ R[z] |
A homogeneous polynomial is a polynomial whose terms do all have the same total degree. For example, h = x5 + 7x2y1z2 + 13x1y2z2 is a homogeneous polynomial of degree 5, in three variables.
The package introduces a concept Polynomial_d, a concept for multivariate polynomials in d variables. Though the concept is written for an arbitrary number of variables, the number of variables is considered as fixed for a particular model of Polynomial_d. The concept also allows univariate polynomials.
First of all a model of Polynomial_d is considered as an algebraic structure, that is, the ring operations {+, -, ⋅ } are provided due to the fact that Polynomial_d refines at least the concept IntegralDomainWithoutDivision. However, a model of Polynomial_d has to be accompanied by a traits class CGAL::Polynomial_traits_d<Polynomial_d> being a model of PolynomialTraits_d. In principal the traits class provides all further functionalities on polynomials.
Given a d-variate polynomial over some base ring R there are at least two different possible views on such a polynomial.
According to these two different views the traits class is required to provide two different coefficient types:
Another important type which is introduced by this package is CGAL::Exponent_vector. It is derived from std::vector<int> and used to identify multivariate monomials. For instance the exponent vector containing the sequence [3,2,4] corresponds to the trivariate monomial x03x12x24. Note that a vector with negative exponents is considered as invalid. However, we decided to in principal allow negative exponents as they may appear as intermediate results, in particular we did not derive from std::vector<unsigned int>.
First of all the concept Polynomial_d requires that the model is constructible from int. This is due to the fact that Polynomial_d refines IntegralDomainWithoutDivision which in turn refines FromIntConstructible. Of course this allows only the construction of constant polynomials.
In general a polynomial is constructed using the functor CGAL::Polynomial_traits_d<Polynomial_d>::Construct_polynomial a model of PolynomialTraits_d::ConstructPolynomial. Basically there are two options:
However, in some cases it might be more convenient to just construct the polynomials representing the different variables and to obtain the final polynomial using algebraic expressions. The most elegant way to construct a certain variable is CGAL::Polynomial_traits_d<Polynomial_d>::Shift being a model of PolynomialTraits_d::Shift.
The following example illustrates different ways to construct a bivariate polynomial:
File: examples/Polynomial/construction.cpp
#include <CGAL/config.h>
#include <CGAL/Polynomial.h>
#include <CGAL/Polynomial_traits_d.h>
#include <CGAL/Polynomial_type_generator.h>
int main(){
CGAL::set_pretty_mode(std::cout);
typedef CGAL::Polynomial_type_generator<int,2>::Type Poly_2;
typedef CGAL::Polynomial_traits_d<Poly_2> PT_2;
typedef PT_2::Coefficient_type Poly_1;
typedef PT_2::Innermost_coefficient_type Integer;
PT_2::Construct_polynomial construct_polynomial;
// constructing a constant polynomial from int
Poly_2 two(2); // = 2
std::cout << "A constant polynomial: " << two << std::endl;
// construction from an iterator range of univariate polynomials
std::list<Poly_1> univariate_coeffs;
univariate_coeffs.push_back(Poly_1(3));
univariate_coeffs.push_back(Poly_1(0));
univariate_coeffs.push_back(Poly_1(5));
Poly_2 F = // 5*y^2 + 3
construct_polynomial(univariate_coeffs.begin(),univariate_coeffs.end());
std::cout << "The bivariate polynomial F: " << F << std::endl;
// construction from an iterator range over monomials
std::list<std::pair<CGAL::Exponent_vector, Integer> > innermost_coeffs;
innermost_coeffs.push_back(std::make_pair(CGAL::Exponent_vector(0,0),-2));
innermost_coeffs.push_back(std::make_pair(CGAL::Exponent_vector(3,5),2));
Poly_2 G = // (2*x^3)*y^5 + (-2)
construct_polynomial(innermost_coeffs.begin(),innermost_coeffs.end());
std::cout << "The bivariate polynomial G: " << G << std::endl;
//construction using shift
PT_2::Shift shift;
Poly_2 x = shift(Poly_2(1),1,0); // 'multiply' 1 by x_0^1
Poly_2 y = shift(Poly_2(1),1,1); // 'multiply' 1 by x_1^1
Poly_2 H = 5 * x * y + 3 * y * y; // = 3*y^2 + (5*x)*y
std::cout << "The bivariate polynomial H: " << H << std::endl;
}
In order to obtain a certain coefficient the traits class provides several functors. Note that the functors do not allow a write access to the coefficients.
File: examples/Polynomial/coefficient_access.cpp
#include <CGAL/config.h>
#include <CGAL/Polynomial.h>
#include <CGAL/Polynomial_traits_d.h>
#include <CGAL/Polynomial_type_generator.h>
int main(){
CGAL::set_pretty_mode(std::cout);
typedef CGAL::Polynomial_type_generator<int,2>::Type Poly_2;
typedef CGAL::Polynomial_traits_d<Poly_2> PT_2;
//construction using shift
Poly_2 x = PT_2::Shift()(Poly_2(1),1,0); // = x^1
Poly_2 y = PT_2::Shift()(Poly_2(1),1,1); // = y^1
Poly_2 F // = (11*x^2 + 5*x)*y^4 + (7*x^2)*y^3
= 11 * CGAL::ipower(y,4) * CGAL::ipower(x,2)
+ 5 * CGAL::ipower(y,4) * CGAL::ipower(x,1)
+ 7 * CGAL::ipower(y,3) * CGAL::ipower(x,2);
std::cout << "The bivariate polynomial F: " << F <<"\n"<< std::endl;
PT_2::Get_coefficient get_coefficient;
std::cout << "Coefficient of y^0: "<< get_coefficient(F,0) << std::endl;
std::cout << "Coefficient of y^1: "<< get_coefficient(F,1) << std::endl;
std::cout << "Coefficient of y^2: "<< get_coefficient(F,2) << std::endl;
std::cout << "Coefficient of y^3: "<< get_coefficient(F,3) << std::endl;
std::cout << "Coefficient of y^4: "<< get_coefficient(F,4) << std::endl;
std::cout << "Coefficient of y^5: "<< get_coefficient(F,5) << std::endl;
std::cout << std::endl;
PT_2::Leading_coefficient lcoeff;
std::cout << "Leading coefficient with respect to y: "
<< lcoeff(F) // = 11*x^2 + 5*x
<< std::endl;
PT_2::Get_innermost_coefficient get_icoeff;
std::cout << "Innermost coefficient of monomial x^1y^4: "
<< get_icoeff(F,CGAL::Exponent_vector(1,4)) // = 5
<< std::endl;
PT_2::Innermost_leading_coefficient ilcoeff;
std::cout << "Innermost leading coefficient with respect to y: "
<< ilcoeff(F) // = 11
<< std::endl;
}
There are three functors in PolynomialTraits_d related to the degree of a polynomial.
File: examples/Polynomial/degree.cpp
#include <CGAL/config.h>
#include <CGAL/Polynomial.h>
#include <CGAL/Polynomial_traits_d.h>
#include <CGAL/Polynomial_type_generator.h>
int main(){
CGAL::set_pretty_mode(std::cout);
typedef CGAL::Polynomial_type_generator<int,2>::Type Poly_2;
typedef CGAL::Polynomial_traits_d<Poly_2> PT_2;
//construction using shift
Poly_2 x = PT_2::Shift()(Poly_2(1),1,0); // x_0^1
Poly_2 y = PT_2::Shift()(Poly_2(1),1,1); // x_1^1
Poly_2 F // = (11*x^2 + 5*x)*y^4 + (7*x^2)*y^3
= 11 * CGAL::ipower(y,4) * CGAL::ipower(x,2)
+ 5 * CGAL::ipower(y,4) * CGAL::ipower(x,1)
+ 7 * CGAL::ipower(y,3) * CGAL::ipower(x,2);
std::cout << "The bivariate polynomial F: " << F <<"\n"<< std::endl;
PT_2::Degree degree;
PT_2::Total_degree total_degree;
PT_2::Degree_vector degree_vector;
std::cout << "The degree of F with respect to y: "<< degree(F) // = 4
<< std::endl;
std::cout << "The degree of F with respect to x: "<< degree(F,0) // = 2
<< std::endl;
std::cout << "The total degree of F : "<< total_degree(F) // = 6
<< std::endl;
std::cout << "The degree vector of F : "<< degree_vector(F)// = (2,4)
<< std::endl;
}
Given for instance a bivariate polynomial it is conceivable that one wants to
interchange the role of x and y. That is one wants to interpret the
x as y and vice versa.
For such a case the polynomial traits provides PolynomialTraits_d::Swap:
Given a polynomial p and to two indices i and j,
the functor returns the polynomial in which xi is substituted by xj and
vice versa, that is, the variables swap their positions.
The order of the other variables remains untouched.
Another scenario is, that a particular variable should be moved to another position, for instance, it should become the outermost variable while the relative order of the other variables remains unchanged. For such a case the polynomial traits provides PolynomialTraits_d::Move.
Of course there is also a general method to interchange the order of variables, namely PolynomialTraits_d::Permute.
File: examples/Polynomial/swap_move.cpp
#include <CGAL/config.h>
#include <CGAL/Polynomial.h>
#include <CGAL/Polynomial_traits_d.h>
#include <CGAL/Polynomial_type_generator.h>
int main(){
CGAL::set_pretty_mode(std::cout);
typedef CGAL::Polynomial_type_generator<int,3>::Type Poly_3;
typedef CGAL::Polynomial_traits_d<Poly_3> PT_3;
//construction using shift
Poly_3 x = PT_3::Shift()(Poly_3(1),1,0); // x_0^1
Poly_3 y = PT_3::Shift()(Poly_3(1),1,1); // x_1^1
Poly_3 z = PT_3::Shift()(Poly_3(1),1,2); // x_2^1
Poly_3 F = x*y*y*z*z*z;
std::cout << "The trivariate polynomial F: " << F << std::endl;
std::cout << std::endl;
PT_3::Swap swap;
PT_3::Move move;
PT_3::Permute permute;
std::cout << "x and z swapped: "<< swap(F,0,2) // = x^3*y^2*z
<< std::endl;
std::cout << "x and y swapped: "<< swap(F,0,1) // = x^2*y*z^3
<< std::endl << std::endl;
std::cout << "x moved to outermost position : "
<< move(F,0,2) // = x^2*y^3*z
<< std::endl;
std::cout << "Same as swap(swap(F,0,1),1,2) : "
<< swap(swap(F,0,1),1,2) // = x^2*y^3*z
<< std::endl;
std::cout << "Same as the permutation (0,1,2)->(2,0,1): ";
std::vector<int> perm;
perm.push_back(2);perm.push_back(0);perm.push_back(1);
std::cout << permute(F,perm.begin(),perm.end())// = x^2*y^3*z
<< std::endl;
}
Since the concept PolynomialTraits_d refines the concept AlgebraicStructureTraits the polynomial traits provides functors for integral division, division with remainder, greatest common divisor, etc. But note that the algebraic structure of a polynomial depends on the algebraic structure of the innermost coefficient, for instance, a gcd is available if and only if the innermost coefficient is a Field or a UniqueFactorizationDomain. Hence, we can not provide a gcd if the innermost coefficient is just an IntegralDomain since it is simply not well defined1. However, if we would consider the polynomial over the quotient field of the integral domain the gcd would be well defined. The only problem is that the result can not be represented over the ring since it contains denominators. Therefore, the PolynomialTraits_d requires functors such as PolynomialTraits_d::GcdUpToConstantFactor. This functor computes the gcd of two polynomials up to a constant factor (utcf). That is, it returns the correct gcd for polynomials over the quotient field, but multiplied by some constant such that the result is representable with coefficients in the ring.
However, note that these 'utcf' functions are usually a bit faster than their strict counterparts. This is due to the fact that the 'utcf' functions are allowed to skip the computation of the correct constant factor. Note that in many cases the constant factor is in fact not needed. In particular if the polynomials are supposed to represent some zero set, that is, an algebraic curve or surface.
The concepts for the related functors are:
Another analog functionality is the pseudo division.
The related functors replace the usual division with remainder in case the
Polynomial is not a EuclideanRing.
The concepts for the related functors are:
File: examples/Polynomial/gcd_up_to_constant_factor.cpp
#include <CGAL/config.h>
#include <CGAL/Polynomial.h>
#include <CGAL/Polynomial_traits_d.h>
#include <CGAL/Polynomial_type_generator.h>
int main(){
CGAL::set_pretty_mode(std::cout);
typedef CGAL::Polynomial_type_generator<int,1>::Type Poly_1;
typedef CGAL::Polynomial_traits_d<Poly_1> PT_1;
PT_1::Shift shift;
PT_1::Gcd gcd;
PT_1::Gcd_up_to_constant_factor gcd_utcf;
PT_1::Multivariate_content mcontent;
PT_1::Canonicalize canonicalize;
//construction using shift
Poly_1 x = shift(Poly_1(1),1,0); // x^1
// common factor 7 * (x^2-2)
Poly_1 F = 21*(x-5)*(x*x-2); // = 21*x^3 + (-105)*x^2 + (-42)*x + 210
Poly_1 G = 14*(x-3)*(x*x-2); // = 14*x^3 + (-42)*x^2 + (-28)*x + 84
std::cout << "The univariate polynomial F: " << F << std::endl;
std::cout << "The univariate polynomial G: " << G << std::endl;
std::cout << "Common multivariate content: "
<< CGAL::gcd(mcontent(F),mcontent(G)) // = 7
<< std::endl;
std::cout << "The gcd of F and G: "
<< gcd(F,G) // = 7*x^2 + (-14)
<< std::endl;
std::cout << "The gcd up to constant factor of F and G: "
<< gcd_utcf(F,G) // = x^2 + (-2)
<< std::endl;
std::cout << "Same as canonicalized gcd of F and G: "
<< canonicalize(gcd_utcf(F,G)) // = x^2 + (-2)
<< std::endl;
}
Of course, it should also be possible to evaluate a polynomial
or substitute its variables. We also require a special functor to
determine whether a polynomial is zero at a given point.
In case the inner most coefficient is RealEmbeddable the traits
also must provide a function to compute the sign at a given point.
The concepts for the related functors are:
The traits is also required to provide variants of these functors that interpret the polynomial as a homogeneous polynomial by adding a virtual homogeneous variable such that each term has the same degree, namely the degree of the polynomial. Of course there is a difference between the univariate and multivariate view. For instance the polynomial
| 5x3 + 7x - 3 |
| 5x3 + 7xw2 -3w3 |
Note that substitute allows the substitution of the variables by any type that is ExplicitInteroperable with the innermost coefficient type. This is a very powerful tool since it allows the substitution of the variables by polynomials. However, for some standard manipulations such as translation or scaling we require special functors since they are expected to be faster than their equivalent implementation using substitution:
File: examples/Polynomial/substitute.cpp
#include <CGAL/config.h>
#include <CGAL/Polynomial.h>
#include <CGAL/Polynomial_traits_d.h>
#include <CGAL/Polynomial_type_generator.h>
int main(){
CGAL::set_pretty_mode(std::cout);
typedef CGAL::Polynomial_type_generator<int,2>::Type Poly_2;
typedef CGAL::Polynomial_traits_d<Poly_2> PT_2;
//construction using shift
Poly_2 x = PT_2::Shift()(Poly_2(1),1,0); // x^1
Poly_2 y = PT_2::Shift()(Poly_2(1),1,1); // y^1
Poly_2 F = 2*x*y + 3*CGAL::ipower(y,3);
std::cout << "The bivariate polynomial F: " << F // = 3*y^3 + (2*x)*y
<< std::endl << std::endl;
PT_2::Evaluate evaluate;
PT_2::Evaluate_homogeneous hevaluate;
// Evaluation considers a polynomials as univariate:
std::cout << "F(5): " << evaluate(F,5) // = 10*x + 375
<< std::endl;
// Evaluate_homogeneous considers F as a homogeneous polynomial in
// the outermost variable only, that is, F is interpreted as
// F(u,v) = 2*x*u*v^2 + 3 * u^3
std::cout << "F(5,7): " << hevaluate(F,5,7) // = 490*x + 375
<< std::endl << std::endl;
PT_2::Substitute substitute;
PT_2::Substitute_homogeneous hsubstitute;
// Substitute considers a polynomials as multivariate, that is, the
// new values for the variables are given by an iterator range
// Note that the value type only has to be interoperable with the innermost
// coefficient
std::list<Poly_2> replacements;
replacements.push_back(x-1); // replace x by x-1
replacements.push_back(y); // replace y by y, i.e., do nothing
std::cout << "The bivariate polynomial F: " << F // = 3*y^3 + (2*x)*y
<< std::endl;
std::cout << "F(x-1,y): " // = 3*y^3 + (2*x + (-2))*y
<< substitute(F,replacements.begin(),replacements.end())
<< std::endl;
// Substitute_homogeneous considers F as a homogeneous polynomial in
// all variable, that is, F is interpreted as
// F(x,y,w) = 2*x*y*w + 3 * y^3
replacements.push_back(y); // replace z by y
std::cout << "F(x-1,y,y): " // = 3*y^3 + (2*x + (-2))*y^2
<< hsubstitute(F,replacements.begin(),replacements.end())
<< std::endl;
}
The PolynomialTraits_d concept also provides more sophisticated functors for computations with polynomials - computing the resultant of two polynomials, their polynomial subresultant sequence, with or without cofactors, and their principal subresultant coefficients.
The principal Sturm-Habicht sequence allows to count the number of real roots of a polynomial using the function
As input, this function requires an iterator range that represents the principal Sturm-Habicht coefficients. This might look complicated at a first sight, as one has to store the principal Sturm-Habicht sequence temporarily. However, we remark an important property of the (principal) Sturm-Habicht sequence. Having a polynomial ft(x) that depends on a parameter t, and its (principal) Sturm-Habicht coefficients stha0(ft), … ,sthan(ft), evaluating stha0(ft) for t=t0 yields a valid (principal) Sturm-Habicht sequence for ft0. The same holds for (principal) subresultants. Thus, it is enough in such situations to compute the sequence once for the parameter t, and call CGAL::number_of_real_roots for each specialized parameter value.We finally remark that computing subresultants and Sturm-Habicht sequences introduces an enormous coefficient blow-up. An application of the functors therefore does not make sense for built-in integers except for toy examples. To avoid overflows, one should use arbitrary size integer types in real applications.
The following example illustrates how two compute resultants of two polynomials, and how to count the number of distinct real roots of a polynomial using its principal Sturm-Habicht coefficients.
File: examples/Polynomial/subresultants.cpp
#include <CGAL/config.h>
#include <CGAL/Polynomial.h>
#include <CGAL/Polynomial_traits_d.h>
#include <CGAL/Polynomial_type_generator.h>
#include <CGAL/Gmpz.h>
int main(){
CGAL::set_pretty_mode(std::cout);
typedef CGAL::Gmpz Int;
typedef CGAL::Polynomial_type_generator<Int,1>::Type Poly_1;
typedef CGAL::Polynomial_traits_d<Poly_1> PT_1;
//construction using shift
Poly_1 x = PT_1::Shift()(Poly_1(1),1); // x^1
Poly_1 F // = (x+1)^2*(x-1)*(2x-1)=2x^4+x^3-3x^2-x+1
= 2 * CGAL::ipower(x,4) + 1 * CGAL::ipower(x,3)
- 3 * CGAL::ipower(x,2) - 1 * CGAL::ipower(x,1)
+ 1 * CGAL::ipower(x,0);
std::cout << "F=" << F << std::endl;
Poly_1 G // = (x+1)*(x+3)=x^2+4*x+3
= 1 * CGAL::ipower(x,2) + 4 * CGAL::ipower(x,1) + 3 * CGAL::ipower(x,0);
std::cout << "G=" << G << std::endl;
// Resultant computation:
PT_1::Resultant resultant;
std::cout << "The resultant of F and G is: " << resultant(F,G) << std::endl;
// It is zero, because F and G have a common factor
// Real root counting:
PT_1::Principal_sturm_habicht_sequence stha;
std::vector<Int> psc;
stha(F,std::back_inserter(psc));
int roots = CGAL::number_of_real_roots(psc.begin(),psc.end());
std::cout << "The number of real roots of F is: " << roots << std::endl; // 3
roots = CGAL::number_of_real_roots(G);
std::cout << "The number of real roots of G is: " << roots << std::endl; // 2
return 0;
}
This package is the result of the integration process of the NumeriX library of Exacus [BEH+05] into Cgal.
The class CGAL::Polynomial<Coeff> had been started by Michael Seel within CGAL as part of the Nef_2 package. As part of the Exacus project it got significantly improved by Arno Eigenwillig and Michael Hemmer.
However, due to the recursive definition the class was rather restricted to the univariate view. Moreover, it is clear that depending on the context other classes that are symmetric in all variables or dedicated for sparse polynomials may be more efficient. As a consequence this package introduced the CGAL::Polynomial_traits_d<Polynomial_d> giving also the symmetric view on polynomials and the opportunity to introduce and use other classes representing polynomials within Cgal.
| 1 | An example for such a number type is the template CGAL::Sqrt_extension<NT,ROOT> representing an algebraic extension of degree two. This is just an IntegralDomain if NT is not a Field. |